library(tidyverse)
library(sensemakr)
library(estimatr)
library(texreg)
library(corrr)
library(lmerTest)
library(lme4)
library(interplot)
library(margins)
library(stargazer)
library(patchwork)
library(countrycode)
library(ggrepel)
library(regweight)

## Functions

zero1 <- function(x, minx=NA, maxx=NA){
  res <- NA
  if(is.na(minx)) res <- (x - min(x,na.rm=T))/(max(x,na.rm=T) -min(x,na.rm=T))
  if(!is.na(minx)) res <- (x - minx)/(maxx -minx)
  res
}

#Load in data

load("cleaned_data/motherdathermerged.RData")

# Study 1 -----------------------------------------------------------------

#Gender Attitude Models

c <- lm_robust(genderatts ~ kinship_score,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")


d <- lm_robust(genderatts ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")

e <- lm_robust(genderatts ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")

texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and Gender Attitudes (Father's Kinship)",scalebox = .9,include.ci=F,file = "output/ess_gender_father.tex",override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = F,use.packages = F)

genderatts <- 
  as.data.frame(rbind(
    summary(c)$coefficients[1,],
    summary(d)$coefficients[1,],
    summary(e)$coefficients[1,]))

## LGBT Attitude Models

c <- lm_robust(gayrights ~ kinship_score,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")


d <- lm_robust(gayrights ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")

e <- lm_robust(gayrights ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")

texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and LGBT Attitudes (Father's Kinship)",scalebox = .9,include.ci=F,file = "output/ess_lgbt_father.tex",override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = FALSE,use.packages = F)

gayatts <- 
  as.data.frame(rbind(
    summary(c)$coefficients[1,],
    summary(d)$coefficients[1,],
    summary(e)$coefficients[1,]))


## Government Benefits Models


c <- lm_robust(benefits ~ kinship_score,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")


d <- lm_robust(benefits ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")

e <- lm_robust(benefits ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")


texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and Government Benefits (Father's Kinship)",scalebox = .9,include.ci=F,file = "output/ess_benfits_father.tex",override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = F,use.packages = F)

benefits <- 
  as.data.frame(rbind(
    summary(c)$coefficients[1,],
    summary(d)$coefficients[1,],
    summary(e)$coefficients[1,]))

## Government Responsibility Models


c <- lm_robust(govtresp ~ kinship_score,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")


d <- lm_robust(govtresp ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")

e <- lm_robust(govtresp ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")

texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and Government Responsibility (Father's Kinship)",scalebox = .9,include.ci=F,file = "output/ess_govtresp_father.tex",override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = F,use.packages = F)

govtresp <- 
  as.data.frame(rbind(
    summary(c)$coefficients[1,],
    summary(d)$coefficients[1,],
    summary(e)$coefficients[1,]))


## Income Differences Models


c <- lm_robust(incomediffs ~ kinship_score,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")


d <- lm_robust(incomediffs ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")

e <- lm_robust(incomediffs ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata), 
               weights=fatherdata$dweight,
               se_type="stata")


texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and Income Differences (Father's Kinship)",scalebox = .9,include.ci=F,file = "output/ess_incomediffs_father.tex",override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = F,use.packages = F)

incomediffs <- 
  as.data.frame(rbind(
    summary(c)$coefficients[1,],
    summary(d)$coefficients[1,],
    summary(e)$coefficients[1,]))


incomediffs$dv="Income Differences"
incomediffs$domain="Economic"
gayatts$dv="Gay Rights"
gayatts$domain="Culture"

govtresp$dv="Government Responsibility"
govtresp$domain="Economic"

genderatts$dv="Women's Role"
genderatts$domain="Culture"

benefits$dv="Government Benefits"
benefits$domain="Economic"

## ESS Main Effects Figure

forplot <- rbind(incomediffs,gayatts,govtresp,genderatts,benefits)
forplot$specification <- rep(1:3,5)

forplot$ymax <- forplot$`CI Upper`
forplot$ymin <- forplot$`CI Lower`

forplot1 <- forplot %>% filter(specification != 4)
a <- ggplot(subset(forplot1,domain=='Culture'),aes(x=as.factor(specification),y=Estimate))+geom_pointrange(aes(ymin=ymin,ymax=ymax))+facet_wrap(~dv)+ylim(-.1,.25)+theme_bw()+geom_hline(yintercept=0,lty=2)+ggtitle("Panel A: Impact of Ancestral Kinship Strength on Cultural Policy Attitudes") + xlab("Specification")

b <- ggplot(subset(forplot1,domain=='Economic'),aes(x=as.factor(specification),y=Estimate))+geom_pointrange(aes(ymin=ymin,ymax=ymax))+facet_wrap(~dv)+ylim(-.1,.25)+theme_bw()+geom_hline(yintercept=0,lty=2)+ggtitle("Panel B: Impact of Ancestral Kinship Strength on Economic Policy Attitudes")+ xlab("Specification")

a/b
ggsave("output/essmaineffects.png",width=8,height=6, dpi = 600, device = 'png')

## Marginal Effects Figures

fatherdata$essround <- as.factor(fatherdata$essround)
a <- lmer(benefits ~ kinship_score*activism+ln_time_obs_ea+small_scale+zero1(agea)+zero1(gndr)+ agriculture + malariaindex+cntry+as.factor(essround)+eisced+(1|country),subset(fatherdata), weights=fatherdata$dweight) 

a1 <- interplot(a,var1 = "kinship_score",var2="activism",hist=T) + labs(caption = "")+xlab("Political Activism")+ggtitle("Panel A. Marginal Effect of Kinship Strength on Government Benefits")+ylab("Coefficient")+theme_bw()
a1
a <- lmer( govtresp ~ kinship_score*activism  +ln_time_obs_ea+small_scale+(agea)+(gndr)+ agriculture + malariaindex+cntry+(essround)+(1|country),subset(fatherdata),weights=fatherdata$dweight)
summary(a)
b1 <- interplot(a,var1 = "kinship_score",var2="activism",hist=T) + labs(caption = "")+xlab("Political Activism")+ggtitle("Panel B. Marginal Effect of Kinship Strength on Government Responsibility")+ylab("Coefficient")+theme_bw()
summary(b1)
a <- lmer( incomediffs ~ kinship_score*activism   +ln_time_obs_ea+small_scale+zero1(agea)+zero1(gndr)+ agriculture + malariaindex+cntry+as.factor(essround)+(1|country),subset(fatherdata), weights=fatherdata$dweight)
summary(a)

c1 <- interplot(a,var1 = "kinship_score",var2="activism",hist = T) + labs(caption = "")+xlab("Political Activism")+ggtitle("Panel C. Marginal Effect of Kinship Strength on Income Differences")+ylab("Coefficient")+theme_bw()

a1/b1/c1
ggsave(file="output/moderationess.png",width=8,height=12, dpi = 600, device = 'png')
ggsave(b1,file="output/moderationess_presentation.png",width=7,height=7, dpi = 600, device = 'png')

af <- lm_robust( benefits ~ kinship_score*activism +ln_time_obs_ea+small_scale+zero1(agea)+zero1(gndr)+ agriculture + malariaindex,data = subset(fatherdata), weights=fatherdata$dweight,fixed_effects  = ~cntry+essround,clusters = country,se_type="CR0") 
summary(af)
bf <- lm_robust( govtresp ~ kinship_score*activism +ln_time_obs_ea+small_scale+zero1(agea)+zero1(gndr)+ agriculture + malariaindex,data = subset(fatherdata), weights=fatherdata$dweight,fixed_effects  = ~cntry+essround,clusters = country,se_type="CR0") 
summary(bf)
cf <- lm_robust(incomediffs ~ kinship_score*activism +ln_time_obs_ea+small_scale+zero1(agea)+zero1(gndr)+ agriculture + malariaindex,data = subset(fatherdata), weights=fatherdata$dweight,fixed_effects  = ~cntry+essround,clusters = country,se_type="CR0") 

summary(cf)
texreg(list(af,bf,cf),custom.coef.names = c("Kinship Score","Political Engagement","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture","Malaria Index","Kinship Score x Activism"),file = "output/interactioness.tex",custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Impact of Kinship Score x Engagement (Father's Kinship)",scalebox = .9,include.ci=F,table = F,use.packages = F,custom.model.names = c("Benefits","Govt Responsibility","Income Diffs"),override.pvalues = list(af$p.value/2,bf$p.value/2,cf$p.value/2))

# Impact of Ancestral Kinship on Closed Personality

c <- lm_robust( nfs ~ kinship_score ,clusters = country,fixed_effects = ~cntry+essround,subset(fatherdata), weights=fatherdata$dweight,se_type="CR0")
d <- lm_robust( nfs ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,clusters = country,fixed_effects = ~cntry+essround,subset(fatherdata), weights=fatherdata$dweight,se_type="CR0")
e <- lm_robust( nfs ~kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,clusters = country,fixed_effects = ~cntry+essround,subset(fatherdata), weights=fatherdata$dweight,se_type="CR0")
summary(e)

texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),file = "output/ess_nfs_father.tex",custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and Needs for Security (Father's Kinship)",scalebox = .9,include.ci=F,override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = F,use.packages = F)
cor(fatherdata$nfs,fatherdata$kinship_score,use="pairwise.complete.obs")

### Sensitivity


fatherdata$islam <- fatherdata$rlgdnm==6

e <- lm(gayrights ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + islam + cntry + as.factor(essround), 
        subset(fatherdata), weights = dweight)
ss <- sensemakr(e,treatment = "kinship_score",benchmark_covariates = "islamTRUE",ky = seq(0,10,by=1))
summary(ss)
plot(ss)

e <- lm(genderatts ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + islam + cntry + as.factor(essround), 
        subset(fatherdata), weights = dweight)
ss <- sensemakr(e,treatment = "kinship_score",benchmark_covariates = "islamTRUE",ky = seq(0,5,by=1))
summary(ss)
plot(ss)

e <- lm(incomediffs ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + islam + cntry + as.factor(essround), 
        subset(fatherdata), weights = dweight)
ss <- sensemakr(e,treatment = "kinship_score",benchmark_covariates = "islamTRUE",ky = seq(0,10,by=1))
summary(ss)
plot(ss)

e <- lm(benefits ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + islam + cntry + as.factor(essround), 
        subset(fatherdata), weights = dweight)
ss <- sensemakr(e,treatment = "kinship_score",benchmark_covariates = "islamTRUE",ky = seq(0,10,by=1))
summary(ss)
plot(ss)

e <- lm(govtresp ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + islam + cntry + as.factor(essround), 
        subset(fatherdata), weights = dweight)
ss <- sensemakr(e,treatment = "kinship_score",benchmark_covariates = "islamTRUE",ky = seq(0,10,by=1))
summary(ss)
plot(ss)

## Mediation Tables
c <- lm_robust( gayrights ~SEcenter +  COcenter + TRcenter + BEcenter + UNcenter + SDcenter + STcenter + HEcenter + ACcenter + POcenter+ kinship_score + ln_time_obs_ea+small_scale+agea+gndr+ cath1900 + muslim1900 + urbfrac1900 + popd1900 + malariaindex ,clusters = country,fixed_effects = ~cntry+essround,subset(fatherdata), weights=fatherdata$dweight,se_type="stata")
summary(c)



# Mother Results ----------------------------------------------------------


motherdata <- motherdata %>% mutate(eisced = as.factor(eisced),
                                    rlgdnm = as.factor(rlgdnm))

## Gender Attitudes

c <- lm_robust(genderatts ~ kinship_score,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")


d <- lm_robust(genderatts ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

e <- lm_robust(genderatts ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

f <- lm_robust(genderatts ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900 + nfs,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and Gender Attitudes (Mother's Kinship)",scalebox = .9,include.ci=F,file = "output/ess_gender_mother.tex",override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = F,use.packages = F)

genderatts <- 
  as.data.frame(rbind(
    summary(c)$coefficients[1,],
    summary(d)$coefficients[1,],
    summary(e)$coefficients[1,],
    summary(f)$coefficients[1,]))


## LGBT Attitudes


c <- lm_robust(gayrights ~ kinship_score,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")


d <- lm_robust(gayrights ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

e <- lm_robust(gayrights ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

f <- lm_robust(gayrights ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900 + nfs,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

stargazer::stargazer(c,d,e,f)

texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and LGBT Attitudes (Mother's Kinship)",scalebox = .9,include.ci=F,file = "output/ess_lgbt_mother.tex",override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = F,use.packages = F)

gayatts <- 
  as.data.frame(rbind(
    summary(c)$coefficients[1,],
    summary(d)$coefficients[1,],
    summary(e)$coefficients[1,],
    summary(f)$coefficients[1,]))


## Government Benefits


c <- lm_robust(benefits ~ kinship_score,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")


d <- lm_robust(benefits ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

e <- lm_robust(benefits ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

f <- lm_robust(benefits ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900 + nfs,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and Government Benefits (Mother's Kinship)",scalebox = .9,include.ci=F,file = "output/ess_benfits_mother.tex",override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = F,use.packages = F)

benefits <- 
  as.data.frame(rbind(
    summary(c)$coefficients[1,],
    summary(d)$coefficients[1,],
    summary(e)$coefficients[1,],
    summary(f)$coefficients[1,]))


## Government Responsibility


c <- lm_robust(govtresp ~ kinship_score,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")


d <- lm_robust(govtresp ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

e <- lm_robust(govtresp ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

f <- lm_robust(govtresp ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900 + nfs,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

stargazer::stargazer(c,d,e,f)

texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and Government Responsibility (Mother's Kinship)",scalebox = .9,include.ci=F,file = "output/ess_govtresp_mother.tex",override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = F,use.packages = F)

govtresp <- 
  as.data.frame(rbind(
    summary(c)$coefficients[1,],
    summary(d)$coefficients[1,],
    summary(e)$coefficients[1,],
    summary(f)$coefficients[1,]))


## Income Differences


c <- lm_robust(incomediffs ~ kinship_score,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")


d <- lm_robust(incomediffs ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

e <- lm_robust(incomediffs ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

f <- lm_robust(incomediffs ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900 + nfs,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(motherdata), 
               weights=motherdata$dweight,
               se_type="stata")

stargazer::stargazer(c,d,e,f)

texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and Income Differences (Mother's Kinship)",scalebox = .9,include.ci=F,file = "output/ess_incomediffs_mother.tex",override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = F,use.packages = F)

incomediffs <- 
  as.data.frame(rbind(
    summary(c)$coefficients[1,],
    summary(d)$coefficients[1,],
    summary(e)$coefficients[1,],
    summary(f)$coefficients[1,]))


incomediffs$dv="Income Differences"
incomediffs$domain="Economic"
gayatts$dv="Gay Rights"
gayatts$domain="Culture"

govtresp$dv="Government Responsibility"
govtresp$domain="Economic"

genderatts$dv="Women's Role"
genderatts$domain="Culture"

benefits$dv="Government Benefits"
benefits$domain="Economic"


forplot <- rbind(incomediffs,gayatts,govtresp,genderatts,benefits)
forplot$specification <- rep(1:4,5)

forplot$ymax <- forplot$`CI Upper`
forplot$ymin <- forplot$`CI Lower`
a <- ggplot(subset(forplot,domain=='Culture'),aes(x=specification,y=Estimate))+geom_pointrange(aes(ymin=ymin,ymax=ymax))+facet_wrap(~dv)+ylim(-.1,.25)+theme_bw()+geom_hline(yintercept=0,lty=2)+ggtitle("Panel A: Impact of Ancestral Kinship Strength on Cultural Policy Attitudes")

b <- ggplot(subset(forplot,domain=='Economic'),aes(x=specification,y=Estimate))+geom_pointrange(aes(ymin=ymin,ymax=ymax))+facet_wrap(~dv)+ylim(-.1,.25)+theme_bw()+geom_hline(yintercept=0,lty=2)+ggtitle("Panel B: Impact of Ancestral Kinship Strength on Economic Policy Attitudes")


a/b
ggsave("output/essmaineffects_mother.png",width=8,height=6, dpi = 600, device = 'png')


motherdata$essround <- as.factor(motherdata$essround)
a <- lmer( govtresp ~ kinship_score*activism +ln_time_obs_ea+small_scale+zero1(agea)+zero1(gndr)+ cntry+as.factor(essround)+eisced+(1|country),subset(motherdata), weights=motherdata$dweight) 


a1 <- interplot(a,var1 = "kinship_score",var2="activism",hist=T) + labs(caption = "")+xlab("Political Activism")+ggtitle("Panel A. Marginal Effect of Kinship Strength on Government Benefits")+ylab("Coefficient")+theme_bw()
a1
a <- lmer( govtresp ~ kinship_score*activism  +ln_time_obs_ea+small_scale+(agea)+(gndr)+ cntry+(essround)+(1|country),subset(motherdata),weights=motherdata$dweight)
summary(a)
b1 <- interplot(a,var1 = "kinship_score",var2="activism",hist=T) + labs(caption = "")+xlab("Political Activism")+ggtitle("Panel B. Marginal Effect of Kinship Strength on Government Responsibility")+ylab("Coefficient")+theme_bw()
summary(b1)
a <- lmer( incomediffs ~ kinship_score*activism   +ln_time_obs_ea+small_scale+zero1(agea)+zero1(gndr)+ cntry+as.factor(essround)+(1|country),subset(motherdata), weights=motherdata$dweight)
summary(a)

c1 <- interplot(a,var1 = "kinship_score",var2="activism",hist = T) + labs(caption = "")+xlab("Political Activism")+ggtitle("Panel C. Marginal Effect of Kinship Strength on Income Differences")+ylab("Coefficient")+theme_bw()

a1/b1/c1
ggsave(file="output/moderationess_mother.png",width=8,height=12, dpi = 600, device = 'png')
ggsave(b1,file="output/moderationess_presentation_mother.png",width=7,height=7, dpi = 600, device = 'png')

af <- lm_robust( benefits ~ kinship_score*activism +ln_time_obs_ea+small_scale+zero1(agea)+zero1(gndr),data = subset(motherdata), weights=motherdata$dweight,fixed_effects  = ~cntry+essround,clusters = country,se_type="CR0") 
summary(af)
bf <- lm_robust( govtresp ~ kinship_score*activism +ln_time_obs_ea+small_scale+zero1(agea)+zero1(gndr),data = subset(motherdata), weights=motherdata$dweight,fixed_effects  = ~cntry+essround,clusters = country,se_type="CR0") 
summary(bf)
cf <- lm_robust(incomediffs ~ kinship_score*activism +ln_time_obs_ea+small_scale+zero1(agea)+zero1(gndr),data = subset(motherdata), weights=motherdata$dweight,fixed_effects  = ~cntry+essround,clusters = country,se_type="CR0") 

summary(cf)
texreg(list(af,bf,cf),custom.coef.names = c("Kinship Score","Political Engagement","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Kinship Score x Activism"),file = "output/interactioness_mother.tex",custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Impact of Kinship Score x Engagement (Mother's Kinship)",scalebox = .9,include.ci=F,table = F)
texreg(list(af,bf,cf),include.ci=F)


# Impact of Ancestral Kinship on Closed Personality  

c <- lm_robust( nfs ~ kinship_score ,clusters = country,fixed_effects = ~cntry+essround,subset(motherdata), weights=motherdata$dweight,se_type="CR0")
d <- lm_robust( nfs ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,clusters = country,fixed_effects = ~cntry+essround,subset(motherdata), weights=motherdata$dweight,se_type="CR0")
e <- lm_robust( nfs ~kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + eisced + rlgdnm + cath1900 + muslim1900 + urbfrac1900 + popd1900,clusters = country,fixed_effects = ~cntry+essround,subset(motherdata), weights=motherdata$dweight,se_type="CR0")
summary(e)

texreg(list(c,d,e),custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Education, less than lower secondary","Education, lower secondary","Education, lower tier upper secondary","Education, upper tier upper secondary","Education, advanced vocational, sub-degree","Education, lower tertiary education, BA level","Education, higher tertiarty education","Education, other","Protestant","Eastern Orthodox","Other Christian","Jewish","Islamic","Eastern religions","Other non-christian","Catholic Population in 1900","Muslim Population in 1900","Fraction living in urban areas in 1900","Population Density in 1900"),file = "output/ess_nfs_mother.tex",custom.note = "All models include country-of-residence  and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship and Needs for Security (Mother's Kinship)",scalebox = .9,include.ci=F,override.pvalues = list(c$p.value/2,d$p.value/2,e$p.value/2),table = F,use.packages = F)
cor(motherdata$nfs,motherdata$kinship_score,use="pairwise.complete.obs")

## Sensitivity

motherdata$islam <- motherdata$rlgdnm==6

e <- lm(gayrights ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + islam + cntry + as.factor(essround), 
        subset(motherdata), weights = dweight)
ss <- sensemakr(e,treatment = "kinship_score",benchmark_covariates = "islamTRUE",ky = seq(0,10,by=1))
summary(ss)
plot(ss)

e <- lm(genderatts ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + islam + cntry + as.factor(essround), 
        subset(motherdata), weights = dweight)
ss <- sensemakr(e,treatment = "kinship_score",benchmark_covariates = "islamTRUE",ky = seq(0,10,by=1))
summary(ss)
plot(ss)

e <- lm(incomediffs ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + islam + cntry + as.factor(essround), 
        subset(motherdata), weights = dweight)
ss <- sensemakr(e,treatment = "kinship_score",benchmark_covariates = "islamTRUE",ky = seq(0,10,by=1))
summary(ss)
plot(ss)

e <- lm(benefits ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + islam + cntry + as.factor(essround), 
        subset(motherdata), weights = dweight)
ss <- sensemakr(e,treatment = "kinship_score",benchmark_covariates = "islamTRUE",ky = seq(0,10,by=1))
summary(ss)
plot(ss)

e <- lm(govtresp ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex + islam + cntry + as.factor(essround), 
        subset(motherdata), weights = dweight)
ss <- sensemakr(e,treatment = "kinship_score",benchmark_covariates = "islamTRUE",ky = seq(0,10,by=1))
summary(ss)
plot(ss)

## Mediation Tables
c <- lm_robust( gayrights ~SEcenter +  COcenter + TRcenter + BEcenter + UNcenter + SDcenter + STcenter + HEcenter + ACcenter + POcenter+ kinship_score + ln_time_obs_ea+small_scale+agea+gndr+ cath1900 + muslim1900 + urbfrac1900 + popd1900 + malariaindex ,clusters = country,fixed_effects = ~cntry+essround,subset(motherdata), weights=motherdata$dweight,se_type="stata")
summary(c)

## ESS

m1a <- lm_robust(opennesstochange ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
                 clusters = country,
                 fixed_effects = ~cntry+essround,subset(fatherdata), 
                 weights=fatherdata$dweight,
                 se_type="stata")
summary(m1a) 

m1b <- lm_robust(conservation ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
                 clusters = country,
                 fixed_effects = ~cntry+essround,subset(fatherdata), 
                 weights=fatherdata$dweight,
                 se_type="stata")
summary(m1b)

m1c <- lm_robust(selfenhancement ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
                 clusters = country,
                 fixed_effects = ~cntry+essround,subset(fatherdata), 
                 weights=fatherdata$dweight,
                 se_type="stata")
summary(m1c)

m1d <- lm_robust(selftranscendence ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
                 clusters = country,
                 fixed_effects = ~cntry+essround,subset(fatherdata), 
                 weights=fatherdata$dweight,
                 se_type="stata")
summary(m1d)


texreg(list(m1a, m1b, m1c, m1d), custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index"), custom.model.names=c("Openness to Change", "Conservation", "Self-Enhancement", "Self-Transcendence"), file = "output/ess_moderation.tex",custom.note = "All models include country-of-residence and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship Predicting Moral Values",scalebox = .9,include.ci=F,override.pvalues = list(m1a$p.value/2,m1b$p.value/2,m1c$p.value/2,m1d$p.value/2),table = F,use.packages = F)


fatherdata1 <- fatherdata %>% drop_na(genderatts, kinship_score, ln_time_obs_ea, small_scale, agea, gndr, agriculture, malariaindex, opennesstochange, conservation, selfenhancement, selftranscendence)

d <- lm_robust(genderatts ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata1), 
               weights=fatherdata1$dweight,
               se_type="stata")

m2a <- lm_robust(genderatts ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex +
                   opennesstochange + conservation + selfenhancement + selftranscendence,
                 clusters = country,
                 fixed_effects = ~cntry+essround,subset(fatherdata1), 
                 weights=fatherdata1$dweight,
                 se_type="stata")
summary(m2a)

fatherdata2 <- fatherdata %>% drop_na(gayrights, kinship_score, ln_time_obs_ea, small_scale, agea, gndr, agriculture, malariaindex, opennesstochange, conservation, selfenhancement, selftranscendence)

e <- lm_robust(gayrights ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
               clusters = country,
               fixed_effects = ~cntry+essround,subset(fatherdata2), 
               weights=fatherdata2$dweight,
               se_type="stata")

m2b <- lm_robust(gayrights ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex +
                   opennesstochange + conservation + selfenhancement + selftranscendence,
                 clusters = country,
                 fixed_effects = ~cntry+essround,subset(fatherdata2), 
                 weights=fatherdata2$dweight,
                 se_type="stata")
summary(m2b)

texreg(list(d,m2a,e,m2b), custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Openness to Change", "Conservation", "Self-Enhancement", "Self-Transcendence"), file = "output/ess_moderation2.tex",custom.note = "All models include country-of-residence and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship Models with Moral Values",scalebox = .9,include.ci=F,override.pvalues = list(d$p.value/2,m2a$p.value/2,e$p.value/2,m2b$p.value/2),table = F,use.packages = F,custom.model.names = c("Gender Attitudes","Gender Attitudes","LGBT Attitudes","LGBT Attitudes"))


fatherdata3 <- fatherdata %>% drop_na(benefits, kinship_score, ln_time_obs_ea, small_scale, agea, gndr, agriculture, malariaindex, opennesstochange, conservation, selfenhancement, selftranscendence)

aa <- lm_robust(benefits ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
                clusters = country,
                fixed_effects = ~cntry+essround,subset(fatherdata3), 
                weights=fatherdata3$dweight,
                se_type="stata")

aa1 <- lm_robust(benefits ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex +
                   opennesstochange + conservation + selfenhancement + selftranscendence,
                 clusters = country,
                 fixed_effects = ~cntry+essround,subset(fatherdata3), 
                 weights=fatherdata3$dweight,
                 se_type="stata")

fatherdata4 <- fatherdata %>% drop_na(govtresp, kinship_score, ln_time_obs_ea, small_scale, agea, gndr, agriculture, malariaindex, opennesstochange, conservation, selfenhancement, selftranscendence)

aa2 <- lm_robust(govtresp ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
                 clusters = country,
                 fixed_effects = ~cntry+essround,subset(fatherdata4), 
                 weights=fatherdata4$dweight,
                 se_type="stata")

aa3 <- lm_robust(govtresp ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex +
                   opennesstochange + conservation + selfenhancement + selftranscendence,
                 clusters = country,
                 fixed_effects = ~cntry+essround,subset(fatherdata4), 
                 weights=fatherdata4$dweight,
                 se_type="stata")

fatherdata5 <- fatherdata %>% drop_na(incomediffs, kinship_score, ln_time_obs_ea, small_scale, agea, gndr, agriculture, malariaindex, opennesstochange, conservation, selfenhancement, selftranscendence)

aa4 <- lm_robust(incomediffs ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex,
                 clusters = country,
                 fixed_effects = ~cntry+essround,subset(fatherdata5), 
                 weights=fatherdata5$dweight,
                 se_type="stata")

aa5 <- lm_robust(incomediffs ~ kinship_score + ln_time_obs_ea + small_scale + agea + gndr + agriculture + malariaindex +
                   opennesstochange + conservation + selfenhancement + selftranscendence,
                 clusters = country,
                 fixed_effects = ~cntry+essround,subset(fatherdata5), 
                 weights=fatherdata5$dweight,
                 se_type="stata")


texreg(list(aa,aa1,aa2,aa3,aa4,aa5), custom.coef.names = c("Kinship Score","Date entered in Ethnographic Atlas","Small Scale Group","Age","Female","Agriculture", "Malaria Index", "Openness to Change", "Conservation", "Self-Enhancement", "Self-Transcendence"), file = "output/ess_moderation3.tex",custom.note = "All models include country-of-residence and wave fixed effects. SEs are clustered by country-of-origin.",caption = "Kinship Models with Moral Values", scalebox = .9, include.ci=F, override.pvalues = list(aa$p.value/2,aa1$p.value/2,aa2$p.value/2,aa3$p.value/2,aa4$p.value/2,aa5$p.value/2),table = F,use.packages = F,custom.model.names = c("Benefits","Benefits","Govt Resp","Govt Resp","Income Diffs","Income Diffs"))

